{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Logistic Regression with Differential Privacy"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We start by importing the required libraries and modules and collecting the data that we need from the [Adult dataset](https://archive.ics.uci.edu/ml/datasets/adult)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import diffprivlib.models as dp\n",
    "import numpy as np\n",
    "from sklearn.linear_model import LogisticRegression"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "X_train = np.loadtxt(\"https://archive.ics.uci.edu/ml/machine-learning-databases/adult/adult.data\",\n",
    "                        usecols=(0, 4, 10, 11, 12), delimiter=\", \")\n",
    "\n",
    "y_train = np.loadtxt(\"https://archive.ics.uci.edu/ml/machine-learning-databases/adult/adult.data\",\n",
    "                        usecols=14, dtype=str, delimiter=\", \")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array(['<=50K', '>50K'], dtype='<U5')"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.unique(y_train)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's also collect the test data from Adult to test our models once they're trained."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "X_test = np.loadtxt(\"https://archive.ics.uci.edu/ml/machine-learning-databases/adult/adult.test\",\n",
    "                        usecols=(0, 4, 10, 11, 12), delimiter=\", \", skiprows=1)\n",
    "\n",
    "y_test = np.loadtxt(\"https://archive.ics.uci.edu/ml/machine-learning-databases/adult/adult.test\",\n",
    "                        usecols=14, dtype=str, delimiter=\", \", skiprows=1)\n",
    "# Must trim trailing period \".\" from label\n",
    "y_test = np.array([a[:-1] for a in y_test])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array(['<=50K', '>50K'], dtype='<U5')"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.unique(y_test)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Logistic Regression with no privacy"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To begin, let's first train a regular (non-private) logistic regression classifier, and test its accuracy."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "LogisticRegression()"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "clf = LogisticRegression(solver=\"lbfgs\")\n",
    "clf.fit(X_train, y_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Non-private test accuracy: 81.25%\n"
     ]
    }
   ],
   "source": [
    "baseline = clf.score(X_test, y_test)\n",
    "print(\"Non-private test accuracy: %.2f%%\" % (baseline * 100))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Differentially private logistic regression"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Using the `diffprivlib.models.LogisticRegression` module of diffprivlib, we can train a logistic regression classifier while satisfying differential privacy.\n",
    "\n",
    "If we don't specify any parameters, the model defaults to `epsilon = 1` and `data_norm = None`. If the norm of the data is not specified at initialisation (as in this case), the norm will be calculated on the data when `.fit()` is first called and a warning will be thrown as it causes a privacy leak. To ensure no additional privacy leakage, we should specify the data norm explicitly as an argument, and choose the bounds indepedently of the data (i.e. using domain knowledge).\n",
    "\n",
    "Additionally, the high `data_norm` that is read from the data in this instance gives poor results, with accuracy only slightly better than random. This is as a result of the large amount of noise requires to protect data spread over a large domain. By clipping the data to a smaller domain, accuracy improves markedly, as demonstrated below."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      ".../site-packages/diffprivlib/models/logistic_regression.py:221: PrivacyLeakWarning: Data norm has not been specified and will be calculated on the data provided.  This will result in additional privacy leakage. To ensure differential privacy and no additional privacy leakage, specify `data_norm` at initialisation.\n",
      "  \"privacy leakage, specify `data_norm` at initialisation.\", PrivacyLeakWarning)\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "LogisticRegression(accountant=BudgetAccountant(spent_budget=[(1.0, 0)]),\n",
       "                   data_norm=99999.04562544584)"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "dp_clf = dp.LogisticRegression()\n",
    "dp_clf.fit(X_train, y_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Differentially private test accuracy (epsilon=1.00): 53.07%\n"
     ]
    }
   ],
   "source": [
    "print(\"Differentially private test accuracy (epsilon=%.2f): %.2f%%\" % \n",
    "     (dp_clf.epsilon, dp_clf.score(X_test, y_test) * 100))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "By setting `epsilon = float(\"inf\")`, we can produce the same result as the non-private logistic regression classifer."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "LogisticRegression(accountant=BudgetAccountant(spent_budget=[(1.0, 0), (inf, 0)]),\n",
       "                   data_norm=100000.0, epsilon=inf)"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "dp_clf = dp.LogisticRegression(epsilon=float(\"inf\"), data_norm=1e5)\n",
    "dp_clf.fit(X_train, y_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Agreement between non-private and differentially private (epsilon=inf) classifiers: 100.00%\n"
     ]
    }
   ],
   "source": [
    "print(\"Agreement between non-private and differentially private (epsilon=inf) classifiers: %.2f%%\" % \n",
    "     (dp_clf.score(X_test, clf.predict(X_test)) * 100))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Tradeoff of accuracy and privacy"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can also visualise the tradeoff between accuracy and `epsilon` using `matplotlib`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "accuracy = []\n",
    "epsilons = np.logspace(-3, 1, 500)\n",
    "\n",
    "for eps in epsilons:\n",
    "    dp_clf = dp.LogisticRegression(epsilon=eps, data_norm=100)\n",
    "    dp_clf.fit(X_train, y_train)\n",
    "    accuracy.append(dp_clf.score(X_test, y_test))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's save the results using `pickle` so we can reproduce the plot easily in the future."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pickle\n",
    "\n",
    "pickle.dump((epsilons, baseline, accuracy), open(\"lr_accuracy_500.p\", \"wb\" ) )"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's plot the results using `matplotlib`. The discontinuty observed near `epsilon = 10` is an artifact of the model. Because of the norm-clipping applied to the dataset before training (`data_norm=100`), the accuracy plateaus without reaching the non-private baseline."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYsAAAEaCAYAAADg2nttAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAgAElEQVR4nOydd5wdVdn4v88t20uS3Wx67wmEkAQSCCWQ0JsgICX60gVfVBT8AaIIFkQpooIFUKIiHV8UaREIBpCYQmgJCYRkU0nv2++95/fHlDt37sy9c3f3srvhfD+f/eydmTMzZ9p5zlPOc0QphUaj0Wg0mQh1dAU0Go1G0/nRwkKj0Wg0WdHCQqPRaDRZ0cJCo9FoNFnRwkKj0Wg0WdHCQqPRaDRZ0cIiICLyOxH5vmP5KhHZLCL7RKRKRKaJyMfm8hc6sq5m/ZaKyPSAZZWIDDd/zxaRH+e1csnzptzTjqK9rllELhSROa3cN/Dz6sqIyHdF5MGOrocmd0SPswARqQV6ATEgDiwD/gzcr5RKeJSPAnuAqUqpd811rwD/UEr98rOqt6M+s4H1SqnvtXJ/BYxQSq1s67E6gna4/jbt39nPp9G0B1qzSHKaUqocGATcDlwP/MGnbC+gCFjqWDfItRwYEYm0Zr+ujIiEO7oO+xvt/R59Ht/L1vK5uFdKqc/9H1ALzHStOxRIAAeYy7OBHwMjgTpAAfuAV4FPzLIN5rpCoBJD2HwKbDD3DZvHugh4E/gFsN3cVgjcCawFNgO/A4rN8tOB9cC1wBbzmBeb264AWoBm89zPuq/JvJa3gF3mvvcCBY5rVcBw53Wavz/AEKJWuSiwDTjY4x5adfyuWaYWuNCxfTbwW+B58/7NdJ3rQ+BUR/kIsBWYaC4/CWwCdgPzgHFZrr8v8LR5jNXANzI8f7se5vLlwEpgB/APoK9j2/HACrMevwH+DVzmeK5vmL/FfL5bMLTQ94EDAj6vsHkfPwH2AouBAR71Hmw+u0sx3pt55vpLzPu5E3gJGJRD/XN5L6uBf2K8VzuA14GQue16jPd+r3m+Geb6W4CHHfU5HaOTtQt4DRjj+i6vA94z6/s4UOTzDIdhfIvbMd6/vwLdHNsHAH8z34ftwL2u5/2hWddlJN85+7vw+DamY7zv12O8l38Bupv3Y6t57/8J9Hfs3wN4CNhobn8m1++sI/+0ZuGDUmoBxstwpGv9R8A4c7GbUupYpdQwjI/pNKVUmVKqCePFigHDgYMxPtLLHIeaAqzC0FJ+gqHNjAQmmPv0A252lO+NIYD6YTQO94lId6XU/Rgfxs/Nc5/mcTlx4FsYH/dhwAzgawFuw5+BWY7lk4FPlVJLfMr3Ns/RD/gf4H4RGeXYfoF5reXAG659HwXOdyyfAGxTSr1tLr8AjABqgLcxrhmv6xeREPAs8K5ZlxnANSJyQrYLFpFjgZ8C5wJ9gDXAY+a2auAp4EagCqMRPNznUMcDR2E800rzeNsDPq9vm/fiZKACo/Gvz1Dto4ExwAkicgaGoDkL6InRgD+aQ/1zeS+vxfhGeprlvwso85lfDRyiDG39BIyGPwURGWnW7RrzGM8Dz4pIgaPYucCJwBBgPIZA80Iwnltf814MwBBMlhb7T4xnOdi8BuuZnmOW+wrGvT4dQ5gEoTeGABiE0QkIYQiDQcBAjM7jvY7yfwFKMNqPGgyhDLl/Zx1DR0urzvCHh2Zhrp8P3KTSexWDMXodEa9jYHw4TZg9MHPd+cBc8/dFwFrHNsHobQ9zrDsMWK2SvZgG1/m2YPhMUuqW7ZrMbdcA/+dY9tMs+mL0tirM5aeA/+dzzOkYwrHUse4J4PuO4/7ZtY/zXMPNc5WYy38FbvY5VzezzpVe14/R4K117XMj8JDP8Zz1+ANGQ25tK8PQBAZjNChvuZ7bOrw1i2OBj4CpmL1tr/P5vD8rgDMCvLeDzfsw1LHuBeBSx3IIQ9AMClj/XN7LHwJ/x9H7djzLLRjaY9S17RZMzQL4PvCEq64bgOmOezLLsf3nwO8CftNfAJY46rwVx/fjKPcS8E2fY2TTLJrx0XTMMhOAnebvPhjWh+4e5QJ/Zx35pzWLzPTDUK9zZRCGKvmpiOwSkV3A7zF6ExbrHL97YvQ4FjvKv2iut9iulIo5lusxGrKsiMhIEfmniGwSkT3AbRgaQEaUUhsxzBJfFJFuwEmYPXofdiql6hzLazA+BIt1+KCUWolhCjhNREoweniPmPUPi8jtIvKJWf9acze/axgE9LXupXk/v4shxLPR16y3Va99GD3Nfua2dY5tCqNn7XU9r2L0Ku8DtojI/SJSEeD8YPSKPwlYFlLv6yDgl47r3oHR6Aetfy7v5R0Y5ro5IrJKRG4wj7sSo0NyC8a1PyYizvfAwn2vE+b5+znKbHL89n3nRaSXeZ4N5jvyMMn3YwCwxvX94NiWy712slUp1eioQ4mI/F5E1ph1mAd0MzWbAcAOpdRO90Fa8Z11CFpY+CAih2C8tG5zSRDWYWgW1UqpbuZfhVJqnKOMcvzehqE5jHOUr1RKBRIGrmN58VtgOUbEUwVGwykBj/0nDBX5HIxe6YYMZbuLSKljeSCGfTZoPS1T1BnAMrPRAcN8dQZGT7USo0cNyWtwH3cdRu+3m+OvXCl1cpbzY9Z3kLVgXk8VRo/3U6C/Y5s4l90opX6llJoEjMUw5XzHp75u1mHY4IPiPN464Kuuay9WSv0nYP0Dv5dKqb1KqWuVUkMxhPu3RWSGue0RpdQRGPdSAT/zqLf7XgtGo5rpHfPjNvM8B5rv+CyS78c6YKCPEzrTva7HEJYWvV3b3c/xWmAUMMWsw1HmekuD62EKAy9y+c46BC0sXIhIhYicimHTfFgp9X6ux1BKfQrMAe4yjxcSkWEicrRP+QTwAPALEakx69EviI3dZDMwNMP2cgwn6z4RGQ1cFfRagGeAicA3MWyr2bhVRApE5EjgVAzHdFAew7D1X4WpVZiUYwjf7Rgf722u/dzXvwDYKyLXi0ixqZkcYHYAsvEocLGITBCRQvNc/1VK1QLPAQeKyBfMhud/SW9AAKOzISJTzDDrOqARwwzhVV83DwI/EpERYjBeRKoC1B0MB/SNIjLOrEelaZcnl/pD9vdSRE4VkeFmI78bwzeWEJFRInKsef8aMQROWgg6hpnyFBGZYd6nazGe838CXquTcoyAgd0i0o+kYAbjffgUuF1ESkWkSESmmdseBK4TkUnmvR4uIpYAewe4wHx/TsTwDWWrQwOwS0R6AD+wNphtwgvAb0Sku4hEReQox765fmefOVpYJHlWRPZi9ABuAu4GLm7D8b4CFGBEV+zEsEP2yVD+egyVfr6pwr6M0UsJwh+Asaap4BmP7ddh9M73Ynz8jwc8LkqpBoyooiEY0SSZ2IRxrRsx1OgrlVLLczjXpxhRW4e76vhnDHPFBoz7Od+1a8r1K6XiGIJqAkYk1DaMRqEyQB1exrClP43RwAwDzjO3bcPo+f0cQ3CNBRZhNHBuKjDu9U6z7tsxzDZp9fXY926MhnQOhpD/A1Ccre5mHf8Poxf/mPkefYBh1si1/haZ3ssR5vI+jOf2G6XUXIwIqtsx7vsmDPPrjR51XYHRm/61WfY0jCCR5iDX6uJWjMZ2N4ZQtN9V8304DcOXshbD9PYlc9uTGI78RzC+j2cwnNZgNNynYURqXWhuy8Q9GM9pG8Y7+qJr+5cx/F/LMXw61zjqmMt31iHoQXmarIjIzcBIpdSsDGWmY2hivmaZ/Q0xoq7WY4QIz+3o+uRKV6///kaQ76wj0ZqFJiOmOn0pcH9H16UzICIniEg308Ri+X7cmk6npavXf3+lK3xneRMWIvJHEdkiIh/4bBcR+ZWIrBSR90RkYr7qomkdInI5hlnuBaXUvI6uTyfhMIzoGcts8gXThNBV6Or13+/oKt9Z3sxQpvNmH0Zs/QEe208Gvo4xAGUK8Eul1JS8VEaj0Wg0bSJvmoUpITONUTgDQ5AopdR8jHjkTA5gjUaj0XQQHemz6EfqAKD1pA7G0Wg0Gk0noUtkShSRKzByr1BaWjpp9OjRHVwjjUaj6VosXrx4m1KqZ/aS3nSksNiAMVrToj8+IzeVkXztfoDJkyerRYsW5b92Go1Gsx8hImuyl/KnI81Q/wC+YkZFTQV2m4OyNBqNRtPJyJtmISKPYmRmrBaR9RhD36MASqnfYaQjPhljdGg9bRstrdFoNJo8kjdhoZQ6P8t2hZGbRqPRaDSdHD2CW6PRaDRZ0cJCo9FoNFnRwkKj0Wg0WdHCQqPRaDRZ0cJCo9FoNFnRwkKj0Wg0WdHCQqPRaDRZ0cJCo9FoNFnRwkKj0Wg0WdHCQqPRaDRZ0cJCo9FoNFnRwkKj0Wg0WdHCQqPRaDRZ0cJCo9FoNFnRwkKj0Wg0WdHCQqPRaDRZ0cJCo9FoNFnRwkKj0Wg0WdHCQqPRaDRZ0cJCo9FoNFnRwkKj0Wg0WdHCQqPRaDRZ0cJCo9FoNFnRwkKj0Wg0WdHCQqPRaDRZ0cJCo9FoNFnRwkKj0Wg0WdHCQqPRaDRZ0cJCo9FoNFnRwkKj0Wg0WdHCQqPRaDRZ0cJCo9FoNFnRwkKj0Wg0WdHCQqPRaDRZiXR0BXJm28fw0CnG74ufS6631un1er1er9fr9d7r20BeNQsROVFEVojIShG5wWP7QBGZKyJLROQ9ETk5n/XRaDQaTesQpVR+DiwSBj4CjgPWAwuB85VSyxxl7geWKKV+KyJjgeeVUoMzHXfy5Mlq0aJFeamzRqPR7K+IyGKl1OTW7p9PzeJQYKVSapVSqhl4DDjDVUYBFebvSmBjHuuj0Wg0mlaST2HRD1jnWF5vrnNyCzBLRNYDzwNf9zqQiFwhIotEZNHWrVvzUVeNRqPRZKCjo6HOB2YrpfoDJwN/EZG0Oiml7ldKTVZKTe7Zs6fvwZZv2sPc5VvyV1tNl+H99bv5aPPejq5Gu5HJXNzYEqe+OZa3cycSiqZYPGOZlniCT7bus+sZTyjiCcWGXQ3srGu2y+XL7N3euOuZSCj7r73Ps31fE/EAx1VKdej9y2c01AZggGO5v7nOyaXAiQBKqbdEpAioBlrV4p94z+sA1N5+SpaSnx3xhOLqR97msiOHMmlQ946uzueG0+59A8j+Liil+NUrKzljQl8GV5dmPe7PXlzOkrU7eeyKw9qlntnqtnTjHnqWFzLltle465yD+OKk/mnlTrhnHmu219vX+vMXl1PfHOeW08cBsKexhRN+MY9bTx/H8eN6+56vsSXOKb96nXMmD2D6qJ70qSymsjjKLc8u5bGF6/jglhPYUdfMzX//gNvOOpDqskL+tWwzu+qbuemZD2iOJXjiq4dxyODunHHfG6zdXs+eRkOI1d5+CkoprvjLYv61bDMHDejGnWeP57f//oTK4ijfOm4khZEQhZEwAAtrd9AST3D4sGr7Xry1ajs/+ueHfHPGcP4yfw0NzXF21DXzuy9PokdJAQtqd7B4zU4OH1bNjX97jxMP6E2/biW8s24nv5s1ib1NMe6bu5IxvSt4aekmbj1jHDXlRXywYTe9K4soL4qws66Fi2cvZOveJiqKI1x6xBAmDerOF+57k3hCEQmFePiyKUwa1J3f/fsT5q/azvZ9zdx17kE8vXg9b6zcxrZ9TZxyYF++e/JoIuEQb3y8jeZ4nLfX7OKV5VtIJBTb65oIh4RYXLG9rpmCcIhZUwexYVc9a7bXM31UDQWREH97ez315nU66VFaQEgEEeheEuXhS6dQURzl8j8vYmh1KQpYvGYnG3c1sLuhhVG9K2gr+XRwRzAc3DMwhMRC4AKl1FJHmReAx5VSs0VkDPAK0E9lqFTFgFHqL8++yhkT3BYtGHyDESI251tHsX5nPU+/vYHn3vs0L8Jjd30Lx/3i39x34UQOGdwjY7mDfjiHG08azVePHtbu9egodtQ189TidVx+5FBEpKOrk4b1LriffSyeYEd9M3OXb+FLhwxk3Y56jvz5XEb2KuPRy6fy+KJ1XHX0MESEpxav5/BhVfTtVpx23CevPIyNuxrs9/DjzXvZsKuBgT1KWLW1jg27Gjjv0AEURsI0xxI8+MYqXv9oG8eOruHyo4b61nvpxt384Y3VrN/ZwILVOwAY2auMjzbvo6a8kHvOm8C4vpV8/dElnDCuF4trd/K3JUYf7L1bjufpxeu59VkjhuTa40aycus+dje08NqKrZwyvg/3XTCRFz/YxKvLN/P1Y0ewdOMexvevpDmWYPqdr6XVp2d5IVv3NgFwzKiezF1hmIFvO/NAzprYj9HffzGl/GkH9eWAvhX89IXlKevPmtiPvpXF3Dt3ZYanBmdP6s+m3Y28sXIbAJMHdacloWhqibN8U+s1xWhYmDyoB2+t2m6vK4yEOOmA3jzzzkYG9CimrslolAvCIZrjibRjnDOpP08uXk9FUYTJg3vwqo8VY/qonry2YiuHDO7Opj2NrNvRkLK9OBrm1PF9CIeEhFIMqS7j4y17+dvbGwgJHDSgG0vW7sp4PSeM60X3kgLeWLmN9TsbePLKw1iwegd3vLQipVy3kii76lsAWPOzU9vk4M6bZqGUionI1cBLQBj4o1JqqYj8EFiklPoHcC3wgIh8C8PZfVEmQQHQFEvw7Sfe9RQWFsf/Yl67XYcfb6zcxpa9TTwwb1VGYdGSMF66WDurrx3N9U+/x7+WbWbiwO5MznD9mXhswVqee/9T/nLplHaunT+TfvwyuxuMj2fiwKSmF08obvjb+/xr2WYOGdyD6rJCrnvyXUoKwiz74Ylpx/nta5/w0ea99nt4nMc79+nuRm44aTR/fquWn79ofMRvrdqeUVjc/sJyXv94G6UFYarLCti2r5mPNu8DYMveJi544L8Mriqhdns98z5K9d+Nv2VOyvJd//ooZbl/92KuengxryzfQnMswXvrd7N8016mDOnBBVMGAlBeGKFbadRu4CxBAaQ01t/9v/d5eP6atPo/++5Gnn03PU7l2Xc30hI3voHff3kSr364hfc27GZkrzIaW+L8d/UOdtW38NTi9QCMqCljUFUpC1Zvt7UTi0MH96BHaQG9K4uoKi3g3fW7qakoRClFTXkRb63azg0njea+V1eyalsdlcVR3lm3yxYUXz16KAtW72DJ2l0srN0JYF/v2D4V/OC0sQzoUcIry7fw5KJ1lBdFOGNCP86dPIDFa3ayalsd8z7aypEjqtm6t8m+L5cfOYSjR9ZwxIhqHnx9FT9+7kN7/e6GFp5/fxMvfPNI+nYrJhxK72CdOr6PIbRH1XDlw4sZUl1KeVGU48b0orggRN9uxRRGwmze02h3YN5cuY0LH/wvSsGyjXtSjvf4FVMZ16+Sh+ev4eiRPRn7s7RT5kReB+UppZ7HcFw7193s+L0MmJbrcTtDP3btjnoABvYoyVguFk/acPcn9pkfcHMsvQcWlBv+9j4Au+qb6VZSkLb9hfc/5anF6/nDRYe0+hyJhOKOOSv48tRB9O1WbAsKgD2NMQrChousMBJmt9kDu+OlFXavvr7Z21b/8Za9tJi9z72NLZ5lXlq6iRtOGs0WR4PrxY66Zs78zZvcevo4Xv94G9ceN5Krjx2OiPDsuxv5+qNLUsrXbq8PcOXp/HfVDt5Zl+yxWo1cQSRk9z5f+850qsoKmXHXa3yytY7qsgK+OKk/Vx8znPKiKCu37GXm3YZgXPbpnvSTeHDu5P6s39nAfz7ZTnlRhBPG9eYElzlsZ10zB//oXwAcOaLa7kA0tsRt7aV3RRFvXH8MkXBmV+u3zP/O9+YX//qIX77yMQDXzBjJzvpmVm+rY9rwar79+Dv8bckGjhhezcOXJTsuX546iC9PHZRy7OtOGMXs/9Ty0EWHUFpoNJ9f+eMCCsLCTaeMtctdesQQ5q7YQllhxF7/wzMOoCga9q33saN72b9nX3yobzmnpmu1hQml2Lyn0V7/pckDmDK0CoAr28mi0fVGcLcCpVS7m0rW7qgDoHdlUcZyVoOyvwmLkPm9xgOaMT/Zuo9LZi/k6asOp7qsEIA+lUV8uruRJWt3cczomrR9rvrr2wDsa4pRVph8VbfsbWTpxj0cMyp9Hzfvrt/Fb1/7hHfW7uLRK6ambGtqiZMw618UDaEwfluCwmLjroaUDxSMnmhVqSHg3vYxGazeVmfXPxNzl29hzfZ6bnve6ImO6l1uv6/u8/pREAllFdxOQTFzTA3rdjSwYvNehvUss4VFZXEUgLIi4/+ImnJuPGmMvV9VaaHnsYujYRpakoJVBKxXIxIOUV5kPL8K87huwuHk9+nsdUcdgmFEr7KsgsKPnuXJehcXhCkuKLbvbX+zwzegR/Z7ffKBfTj5wD4p6/50cXpnRkR42KUxZxIUrca8VUrB5r1JYdGtxPs+t4WOjoZqFaEcG35L/QXYsqeROsfH+4O/f8CTi9Z57ZaRNWbvLpGlsbTMT/udsDCfQdDLemDeKtZsr+dfyzbb68b0MZxuS9bu9Nynh9kY15qNrsURP5vLxQ8t9IwMaWyJ883Hkj3xppi/sG6Mxe2GvCgaxu9R7qxv9lxv2bW37/PXHJpjCVsLs9i8p5Gpt72SUk9IvqfRSPKz7Nstc2fEoro0XTPLRGVxAQ+ZjVxzPMGuhmbKCyN2Y1xWaDRsZUWp/cmSQu8Gr9xVzhI6AJGQ2ELCXc4i7PimI6FUwWFtiniYboJidVC8KDYbccu5nisi4tkZ9VvfnljfoVKKLXuS72FxQfsLpi4pLHK1QznD/g697RXO/M2b9vKf3lrDd556L+cqbN9nNCBOQeRFLN51fBYvLd2UYgffuKuBT3c3eJYVW1jkdl3ORxc1e5OPL1pn96qd9O9u9PRqtyeFxb6mmN2D9rr3//5oK39/J2kztzS7SDj9pWlqSVDfZLwbhZEQflfS4GOKso5tRap4NWbb65rSNIsVm/ayaU8jf39nY4rAs66r0NF7rikPJiyqMjSGXlQUR+jbrZi+lUU88t+1PPRmLZWO3qilyZUXpjbufg2qW6ikCosQ5aaw8NUsHPfO3RmMmmpsONT65qpnub8wPWhAJQBHDK9u9fE7CutO7WpooSmW4JqZIzj/0IFcPG1Iu5+rSwqLXDsYTS713HIYtoW4CqYxtNg+i9bb9j8rvvqXxXzljwvs5cNvf5XDfvqqZ9mQrf4GExZWMWc7YPlzNu9p4v55q9L26VtpCguHZrF6a/K3V+z/LpcWYDXoUQ/zRVMsYWuZhmbhfS1+fgur/jvrmwmHxNPvsnVvurDYXpfsAbbElS2kLE3FqVl4OUK9qCrLTbOwGu8Cx7mcDbxlj3cLAd/jFfoLi2hYqCg2tvtpJk5B6xbs1j2Iegj8oPQs8xe6hw+rZsF3ZzBzbC/fMp2VkHlvNu02TFDDepbx07MOTLn/7Xaudj/iZ4DkqFq4hUV7YGsMHiF2KeUSlhmk3avQodhmqIDXZTWJzmfnFZ7oxGo0Pt6SFO7N8WTD7WWj31mf6my2TTseDc2+ppjtpC2MhHxNak5h4TxOLGEMktpR10L3kgLPc2zb15RmhrK0UkjVWizBVuASbF89amjWXm+PHM1QFaYQcApR53diaRZlhQGFhUtjcDZWYYcZyu/LDaeYnlKv33oPWuuvAKjOoFkA1FQE0+A6G9Zd22aaQnPtNOTC58LB3dSSefRpa7BNIfuJZuEX0eNHKEczlH2bzLf77+9syDrC2mo8F6zeYQcpNLUk76NXJ2CXS1hYvf+Ihwnjzjkr7PKFkbC/Gaol2dgXR8O0xJPLLXHFzrpmupdEPTWQrXub2OO6t9ucwqIlXfg5e/sAN548hvmrtttjD7zI1d5e7iEsnCOtrfWlAYWFW6g4zU1OB7cfIkJIjPfEbc6z6hJtg8+ipMA4/5QhrQvz7qxY5mDbhBnJX/+/SwqLtpqh2oPmgFFOXcVn4Qy7C4L1DIJeli1TlGG6+uZj72TdxxK0n+5uZO2OegZVlaY8Sy/Nwm2GskxVXj4Lp2AJhZyVTKWh2TjPayu2pJmU/v7OBnbUNdO9tIDmePo99DRDORziDY6ILOt6vExmxVkiaXI10XiZoXbUO9NyGP/dWo4fbnOVs9GKhsQWFpkcvpGQMRjO7bOwhIfXM8yFt248Ni/mmY7EulVWx6otfp1sdE0zVI4RBq0VFvXNMVZu8fZvNNlO1mxmqK4RDfXp7lyFRTIKw8mryzcz+IbnWL8zdSyAZYZqiidobPG+Z+68Oy3xhP0xrNyyj9Xb6nh0wVp7u9dzdUcuWT33bI1eLK4ymKFiJBKKix5aaJexnO/feeo9FtTuoEdJQZp/obwwwrZ9zWkO8u11qWaoFleHwqt3mC3sMqhvw8JqgJ335bxDBuZ0DCcl7ugbR3XCYQnUiFlF/DSLtpihAPpUFtsaxv6CdacsC0dbIsay0UWFRW7lLTNUrg32V/+ymJl3/9szeZjVC4xliYZyNwSdlZyFhfnmuC/ryUXGCFxnTD+AZeNZXLsjJbrJSYvLVNcST9g9waZYgtN//QZzHKG3Xg5ut8/Caqi9eutO4gnla1JraI6n+FfOPLgfV7hGYXcvLUj7UCuKo+w2o1ScpGoWsTQNyW2GAmMcSCai4VBG/4K7btaVWuc6Yng1P/7CAY7tub2vbm3A6ZuKBuztWqbCsEuDsDSKtpih9les+x6zNQstLFLI9XZYH2ssR7/B6x8bNmJ3I6KUshuPbELAEibtna2yvbHCP70aKjBexrWOkcOWducelGdFZ7gFs7X0zDsbOemXr3uewx0K2xJXlJo9waZYnL0uc46XGcrtTG5s8TdDpZwroXw10PqWeMq20sJwWu+/vCiS1nsuKQjbjkcnO+qb6WcOCGtoTqRdh5dgC6JZLLxpJkN9kiFaz3VwlTEAbaw5xsUyX3UvTdWMrHEPoXZofMIhYcrQHgzrWcq3jxuZsZzz3BZJM1SXbK7yStIMpTWLdqHJRwsI2oC7G8SWuLJtukGjoYJoFnsbW7j31Y87xGRl9cD9XrU75qzgqDvmsnGXMe4i7OrRWIRt81Tq/kGuyX2slnjCtnU3eZiuvBp3t1nQ6XTOVId4IuE7nqKhOZ7SoEdCoTSHeTQsaR9qqUO49CYAACAASURBVGmGSqt3S8IeYdvgEkTgo1lkcWBHQ2KOTPYuZx3z9IP6Unv7KQwwRy1bgsltprv62OGcM6k/5x0ygNbg1EyiYSMa6pVrp3NAv0rffWxh4WeG0ppFGpZm0aI1C29y7e1YvUt3gx3UNORWSJwmiWzHaMkhN9RPX1jOnXM+Ys7STYHq1Z5Y96glnkjzQyileOsTIwmblefIegRuARz20SyyzYcA6aG0LXFlm1YaPSLavDQL9/OwhEUsrnz9S8N6lhKLq5TIpNRjxFLqFg1LylgIY10o5UO97IghlBVGPDWLeELZQrChJZ5uhvLoQRdmMUNl8glMG17FoWayx0KXhmIJEffxu5UUcMc5B3lGQz108SFcfmTwQV9BNQLr/rmFQjJ0VgsLP+zBp9rBnUqrzVCuxiJoD96tWTg/7qxmqERmE9g763bZEzZZoYsdYbGyGuOEMhrpmEsgJoWAsd7u0biuy1KL3aa7Bh+ntpNYXLFyyz7bJNYST9hRNl5ahJcAcje8VthqLOEtLKJho0ceTyh/zaIlkRJ+HQ2H0uzn0XDIbuRmjqnhe6eOpbQwnJYKJBwSU2MyNItGlz/EqpObbCGRVkPq5c/75XkH271z93EKfNZn4phRNSlJ87w4+YBk/qSgvd2Ij2ZhCcJ8NoRdlaTPwvje3P6edj1X3o6cR3LNDWU1Km7hENSH4d4vRVhkMUMlNQvv7V+4700unr0w5bjOxmJfU8wendkerN1e7zlS2dmrboqlNmBNsYT9IVsvpThe0nhCMfnHL/P04vW2GcotLLw0Azct8QQz7/43J5s+jZZ4wtYsvIWFl2aRum6v6cOIJRKewQiRUIiwGbLpp1k0uDSLSDiU5leIhiXNjFJaEEkT/OGQEHNoFvXNqQ7ugnDIN89QJjKZaAoiIft5uDWLpBBpv1xC910wMWU0dNCwXuu7dgsLa6ktI7j3V6zXotnWLLSwSMHru5l59799y1v2bvcAuqCahdu34fy4s+eGSh2Ul0gobnj6vfRoIbzTPRzxs1eZ+tNXAtUzG0s37uaoO+by0Ju1adsaXIPdnNfY1BK3e3WWJmW1lS3xBI0tcbbta+KmZ963P/Rm130JMjDSupebzDEfLbEERdEwkZDQFIunNSLePovU81qDDf3MUJGQEA1JSnJJJwWREPXN8RSfSYGHGSoSCtm9e1tYeJhwIubsaNagtYaWRIqG1NoG0TL1eGU3KAiHbB9Skdt8FjFDaNtxMJf7WEE1Atvc5NPgaQd3OmmahRYW2fEbDwGOzKNxt2bh39D/5Lll9u80M1Tc6TTNzcG9dV8Tjy1cxxV/XpRW1m7MlNHIbdzVkDYiuS2sM+fgcKfghtSef6PLjt4US9gfsiXQ7Jc0oexrS6ikP8ktHPzGVjhxm4FaEopoOERhJERTSyItSiaIg3tPQ8xe7zXaPmJqBO6QW4uKogj1LlNRxMsMFQnZ5hLrv5ewCIeElkSCUjNHktvB3dpGO6NmEU6mX3drFhZtHfk7vn/Sce32LQTt7YayRGBpB3c61idhtTNas3AhIvz9nQ1245cNq+fmtq9n0iweeH21/dutWTgbvlwd3JbD0yvfvFW22TTHHH67dxK/uqYYP3luGW+u3Jajico/RUdjihkqkTZSOuwSAuKIhnIKOetddZudGh2950uP8HaO7mpITwIYDQuF0TBNsQTuDmqqhpfgln8stc1OFnubDCEQTyhaPIRL2NQI/NKMdy8pYHtdU4pmEQl5OLhDgtVGWh9smUfSPMGIFIuGQxRFQ2mCubXCIlOPMhQS2xxW4GrIrVehLcLib187nDMPTs5c6XbQ56oR+DV42cbKfB6xB+VpzcIbAb752Duc9dv/BCrvN6dB0Ggop2bR2BLnoocWJo+R46A8K+d8d1eG0kRCOVJvJ9i8x3+OhPvnreKB11dz4YP/5eg75ga6Bkj2QrxqnCIsWlKFhdNnceXDbzN3+RZbKBjOcEuzSIYUWwJ1V30zP33+w5TxD5cdOYTTD+qbVoe0JICxRFKziMXTfFVO881rK7Yy+z+1ace0Zr9rSSh7pjQnUXN0sTV15/8cNoifnJkcnDZteDXrdjSwfFNyVriCiJfPIqlZWPX0Gi1svQfRcIjiaJj65ljKvW5tg5jNfJXM+ptaLtEOwmK0Y7Im8IhmyrEBc0d2We+rjoZKx+602ZqFjobyZGuW6SotnGGhTtxmKT/iZkP+3f97nycXr08Jh8zmJLcc4JZ2YuVgcguLhpZk2odsKUScKa5zSWXil6LDOr/VK2mMuc1Q8ZSX8JZnl9qNTCyR1CyUo+6Ws/h7z3zA7+etSklxUVoY8WxA0tKLO81QsXQzlLOO9c3ePgdLCHz46R7+b8mGtO1h02dh8cVJ/TlnUnJswYkHGNN/vuQIZzbGWaSHd9qDx2zNIl1Y2IMEQ0JJgWniCqhZ9M6QGdVqYP384NYzdwtca31bcgpl0yRybeR9fRbaDJWGPSgvpjULT9wKQbY5FT7dZTTQ2aKhFqze4TkJzwcb9nD+A/N55L9r+f4zH9jrSwvCWR3c1vaYLSy8zVANLUm7eKbpMf/+zgYenr82ZZ0x4c8GvvbXxRnrYr1GXrersSVBNyu1RksiLRrKGZJXHA3bwi8WT45mVyqpHS1YvYORN73A/FXb085VWhDxjO7ZWZfULJQyHNIFYaEwEqapJZFiyy4Ih1IEZbZoKz+To3t8RLfi1LQdA3uUEAlJSjqUSFjSGvVoOGTfI+u/02fxwjePZHTvcvvdjYRD1FQUsnFXQ1o0lB/zvzuDYT29R2hbmoVfU2GZHt3bveYZyRW3cEgTHjkKIrfPwlrSobPpuEPYtc/ChVs4ZLMmrdpmOL+dDfvKLXvTTD3n/v4tz0l4rnx4MYvXpE/9WVwQyZ51NpFqAnPOk+vEmVDOHUnkxCtb6//8cQHzV+3gVXO8xlf+uIBfe5hcko5oD82iOW7PlNbk1ixaEikvYWE0bJvmWhxRRorkPX5/w26a4wnPEczhkNhmLGfDstbhg4onDJNWJByiMBpKi4YqjIT4xzsbOes3b9LYEs/qQPcTJuGQpPR8K4ujKY1VJCQURcPsbkgKsgLP0NmktmFpQM6JfobXlKXMMx4NC8N6lrFqa12KYM6HzwKSphy3ULAc3+3ZxqQ5uNtLs9BmqDSSPgsj6WZ7pGfxo0umYHQ3dZnmVBhSXUrt9nriCWU32CGBmXfP890naBqQ0sJw9qyzLs3CMp25B2I5R/J6OWLjjoFxXuxritkN9byPtjLvo618fcaIlDLNjkbdTVMsTq+KQvN3gkgo1QzlNAG9u24X75qhv7FEwlaBlQpuFrMOFw0LVhDU02+vt7cnJy1KmqGcl18QCbFhVwMbdjXw8Pw1Wf1PfrPdRUKSYsJxz7sQMc/vNKM5TU7Ode5xFgPMzLTWdaTMBhcKMbRnKU8tXp+SsTVbdly/q7R63cYkPrvTto/rW8FrK7bSy2XKsm5brhOKATx55WHM/yRdc3T7T3L3WWgHd1CcobP5NtN1ybufrln4NxSje5fTHEuwYWeD7T/INqjPHSrrR0lBJHA0lCWALNu623xlaBZWb907JNRt0wc45UBjpOzexhYj4ieD8LKEkFeVG5rj9rSgjS3xFOdxc8w77BQMIegUfC9/uNmzHMD0UT353iljgOQz8IuUqTPvU9Q0QzW2xFPq7XTI/nf1Dl9hkI1IOOmzqCiKpvXMwqZmkbJPKJSmARQ4hIH1f3hNedqxnOcdWl0GQK0jQWNWzcLndbPOecfZ4zlkcPe07d+aOZJn/ndaWm4m+1VvRTtzyOAeaR0Soy7pWlcu+AkL7bNIxzmfRT79FdBVhYV7OUN7PWmQ8eE8+MYqu2HPpqoFHaxXaqaJyIQ73YdlLmmOxVOEXn1z3G7MvXwWO+qa7elF//A/k+31fSqNnqKVKmTETS/41sXpW7DYuKuBxxeupTGW9Fl4jbPwG6meEjqbhRPG9eayI43U3tZL7vcorPElBZGkZuGsk3O8QFMswZ6G1o1HCYeSPguviXEiIUnLmxQOSVoD6NRQnB/tS9ccxa/PP9gu4yw/sldZ2vla23u2TDTdSgq4YMpAj+0hJgzolrZe+fgy2oI7rDjXRkxrFrnTEldpASDtTZe8+27hkEmzGNunglPH9+G59z61G+xsNzVoSG1JYSRrQ+keZ5GMzErtkW/Y1WBHEHnNTX347a/yL3MuB6cpobrcMB05zSReJBLK7n07b9c3Hl3C9U+/Tzyh6N+jhHBIWL+zIaUO8z7amuJPcJIpQZ8b55wMloPb71bvNsdcREKWzyJVWIzqley1N7bEMwqLUp9MrBZWQ+spLEzNxklCqTT7eTTi8Fk4GrtRvcs5zQwTdkYcRcIhBlclndXWZEpLNyZDdL3IZoaCVM35r5dNCXS8XFPoZCI9b5b2WeQLq+Mb05qFN07h8OGne5hw67/Sykwf1ZNrjxvJpMHdqSotIK5U4CHxmUJqTx3fx7YrF0dDaeMslm3cw9KNSZuxe1rVBkcYr9O+/5PnltkPvq7J26Qy76OtgBGSeaHZe7Sm29zu4Uh2OnVvfXYpP/qnMSrdmT56h0PIVBZHGdijhJVb9vHGx8n5nv+2ZINvI9aSxfTlxDk1qPUE/AS9pVk4zVBOAXbxtMH276aWeNo81066lWSYxF4lfUFe80RHQqG0iYcSCZXmW4g6NJQgoZ/RkKRouD8/ezwAhw+r8q9rBpwNqSU4epYXMm14dcb9rpk5gsmDuqfkcmorbtNi0LDcZBiv9/3Ld2PYFUk6uFXe06F0SWHhbF9qt9V59sS7FUf5+owRFEbCiAgJR1oKr06U0zTjdTyLey+YyM/OPpBuJVEqi6NpWsjJv3qdU371hr3snla1qSVpanKOCt5Z38Lo3uWUFUZSIm+cWKOTSwsj/OTMA6m9/RTbxu2VBM/SNu6f9wl/emuNvb6+OemTcIZ3FkVDDOtZxisfbuHJxYaz+f++drjvvQBDGDbHgmliTtNRcsxHapnqMqNht+6BZYZav9OYR+O4sb3493emM2VoFYu/N5MTxvWisSXhe88Aupf6z7usSDaufuk53HNJJJTH1J8RSZrWAjR21od91kRj5PPkQT14/5bjueOc8b51Bf8wcWd9cmkzBlWV8tRVh7fr3NStTfdhoc1QwXGGzmrNwgPnB+NnMnKq1WEz3UEyCV76TXVnXc3EmQf3552bj6coGs46KM/qdaeboRJp5yktiBANi2/DZ5lanIO9MkXPWCksbnt+ecr6JWt32QLNStonAiNqyhlWU2oLy8riKGP6VNhOdO/ry8EMFXEKC+O/u/G774KJQHI0dzQcYvmmvfb2Qwf3YJBpvqkqK6QoGqYxFs+YQ8s9ANKJcjT8ftOSujWLuEqPTHOagYKYUazft581nsXfm0lBJER5UTRr9ldfM5TjPWjLALv2IE3rCtjIW6ZJ9/1L+re0ZuHGzsrg0YFpb7p86Kyfg9k56CskhrnDnqfW46Ube/NL9u+g4Z+RULoZyk3M5bNocPgsnnvv05SyZUURouGQr/19b1PMjM5JfnyZJsXZurfJtye6css+GlvibN3bxLdmjuSq6cPMNBbCuh31/PSs8VQUGYPnMp3DOYI7G14+C3fkmRV5ZGWLjYSECQO62Vl63dFCRaaJKpOwyGSGUih7EJ2fsHA34EqptEGF0bDY4ad+jXU4xQxlzlAXCVFVVuhbv6C0VrPIB+5GK2iPtz1Gk3/ecL6GWrPwwNm++GsWjt8hIZ5ihsp8U4PMvQBGA+HULP5t+hScOLPOKqXsY+9pbOGnL6T2+MsLIxREQmkJ9ZyUFoRT6p9Js9i4u5F9Pqm3ATsJYf/uxXYjfPDA7vzmwklUFkft81i+hokDu6U4lsFIahhcWKT3mt2yzCpj5ZKKhkPcdMoYDh1izPSWJiyiIfY0xDKaoXpmaIydPTIvM5R1DuvcJx/Ym+PHGilAXv720XYZZ+85iM+itc5av1gO5/E6urFNnxY1VzNUe9Zm/8Y5PkaPs/DA6RT1SxHufGFDIiiVPemfRVDNIhwS+5gbdjXwP39ckFbGGQ3VHE/Y0T9WjignZUURCsKhjA2fu/frFZf/v8cMIxoWNuxsSHFgu7HO45UB14klLKrLCqmpSG14F9bu5N65KzPub1Hk5bNIK2Nczy6zbkXRMNGw4UuB9DEoRdGw76RFANefODqlsXInMFQqGZHllSXWWe+h1aX85sJJ9jzXw2uSoa/OxtqvLxIOYKpqLSmaRQeZa8b1rQDSO2O5h852yWapQ9CaRRacDYyfZuF8YcMixJXKOveERZNH6ogvTx3EU1celrIuEg7ZGoPXBEnH3vWarW3EEyolJYVXeoqyQsMMlUlYuHu/XsLisiOG0qeymN/9+xP+bDq2TxiXHu1iRRBVZHFuWo2jn+153Y4G331/etaB9m+nGWriICPmf3TvVE3FapgtIWcJR8vx7Y768pufAeDlbx/FVdOHOWauCzNjTE1KGUVSk/TSfJzre5T6m7OMGe58NwNuzaJ1n57y9Vok6aiooUcun8oL3zwybX2ujmk9+C44Tj9OvnNndUlh4fxe/HwWKWYo02eRLemfhZeDe+bYXkw2J723sOLJveqweU8jq7bWpdQz22xxpYURohHJmOeozBXe6eUQjUZCdiP7hzdWA3CERwilFV1lzdrmh9VYKlTGAZBu+lQWcf6hyQFiTgf3qeP78sb1x/CXS6bwhQl908pYkVxWfqUBPUqAdCen2/nsVW+rYS6KhtMaUsM0mEgp78YaLZ5JWARp/FPDW9vXDOWko4SFFRDhJtdrdT/jGnMsUVsnaNofcd6pfD/3LurgTobA+pmWnC9cKGSaoWzNIvMXd9XDb6etcw80gmR20Zgj75QfsURyjueSgrBneopy0wyVCbcZyusDKgiHWLk1debA/t1L0spZjvSK4syvgdWIxuIqUM/Wwq31FbsGx1l1uue8g3nmnY1A0mG/w0zDXmrOCfHFif2pb4pxnkP4AGlhrU6snpb17MIhSTPRKJWc0MlXWATQLJymLl+/gtPBnUfDvNVodHT/fGyfCpZ9uifnRsztz/n52QcxY/SmtFQlGpdmkedBi3kV1SJyooisEJGVInKDT5lzRWSZiCwVkUeCHDeZ/CyoZmEsWCOAs/XOvGzgXj1HK6KlJZ7IOuo7oZLHdffkrURylhkqE6UF2c1Q0bBw1zkHMbQ6fYSwE8sMVZ5Fs7B8FrGEwmnJmza8im965Aa6zJwJz22ay9YznD6qp11mx75UzSIcEi6aNiStQbeWvUxA9iA5854aKTlcwgJlz+Lnp6VkC60F413I9qmmjuDO34fdWQav/fWyKTx55WFZA0rcuOtfWRzl3EMG+JT+nLM/+CxEJAzcB5wEjAXOF5GxrjIjgBuBaUqpccA1QY6dCDLOIiWU0CUsgl6EA68HYX3wcQ/NwmtksjUngrMn/8jlUzh4oGG/LzWjoTLh9ll4NcAiwmkH9eXV66ZzybQh9KoopNoVEVRWGGFPQ4yQZE+HUVyQFIpOzSKRMOZ7cHLfBRM5fpwRLeR+NpkajeU/OpEHvzLZCNWNhKgzNa+SDD4JSDbwVaWFaXmmrGdm9fpDXsJCJf1HflqKldTQL1rKOnY2Ukdwt9Jn4Xqt7v/yJK44aqidlNB9no6ke2kBh7hMt0HoKAd9V8T5qLtyNNShwEql1CqlVDPwGHCGq8zlwH1KqZ0ASqktQQ6sUnwW3vZ9p3pm/WzySKQXFK/wP+vhtMRV2uA8L//IxeZ0rFZPvqIowuHDqm3TU0jSE9S5caekyCZcbj5tLPNvnJEytwIYgm5vYwsVjhBZPyzNwppjwkKh6FmeKoQiYbE1paAJGcHQEJy+BTAEYTZfgFW2prwwzdZtfTtDTA1r/c6GtEY9iIO7vilpPmwL7qyz7cHQnmV89+QxKdfV1QevdRbNqCuQEsjThYVFP2CdY3m9uc7JSGCkiLwpIvNF5MRcTxJonIXLDBU0hNaJt2YRMuuQyDoLn5MKs8G3opAsJ3UsnnAIDu99+3ZLnY8gm48DjBfKWe6yI4bQHEuwpzGW1bkNSZt9LJHqsUgo0jSWYT1L7R54ttHtfljaQqaevLtsTUVSWAw1Z5OzntlRI3va5dN9FoojRhjO/0FVSS3pNxdO5MdfMObiDqJZBMGddbY1uDs6mTTerkpXr/9nSapmkd8AgI52cEeAEcB0oD8wT0QOVErtchYSkSuAKwAKeg93rs/gs0iPO3/ozVogc+4nP7x6/NaHGourNKGVSSBZmoXlMP3B6WMpjIY4ZnQNL3xgzPVcXVbIFo85xgf2SJ1WM+jMaikD+SIhmmMJdtQ1eybPc2P7LOKJlBdSKUV1edLpu+CmGdSUF7HFHENiyYqhPUtTIsOyYfXwg/TkLdNRz7JC23770EWHsG5Hg32fo+EQT3z1MKJhSQssUMCVRw3jnEkDUrSkkx0pTo4a0ZO/vb3BM8V3LnjlhmorXjKnq2sWncWM1hVwDsrr8GgoETkNeE4plWsLuwFweqX6m+ucrAf+q5RqAVaLyEcYwmOhs5BS6n7gfoDCPiNSWuEgPgv3t5OLecTC6wWOZoiGyjSq2Wqg+3UznM59Kov55XnGnAeWycSvoXT7CCytpLwwwt4Mo7WdFESM8SH//mirXYdMFDnMUCFx+CwU9HCk0qgpN7SeEpdm8fw3jgw80BGSfhi3M9+zrGWGqkj6LHqUFtj5oyysEeBvmTO7GeHUgDLeFbc5zckXDu7HjDE1WQMBLPwixrxyQ+VKMqW4f3r3rt7YdnVh91kiKZpFx5uhvgR8LCI/F5HRORx7ITBCRIaISAFwHvAPV5lnMLQKRKQawyyVPgl2BoJEQ7klbtD5Kpx4qXjWuj+/Vcs2hxbQvSQ9G60Ta6xEn8r0hnqC6ey2nOFuBlalCgtLYHXPENbpxqmNHBcgNbXTN+O8KqW80yIn04MYE08VRcM5ZTW1NQufEdWpZU0zVHmR3fPP5PexGmm/EeR+BBEU2dq4qtKkQGqrgztTg9rVbf65Rk99nnHeqnzOvw0BNAul1CwRqQDOB2aLiAIeAh5VSu3NsF9MRK4GXgLCwB+VUktF5IfAIqXUP8xtx4vIMiAOfEcplT6pbwYCjbNoh5fPqydoNaIPvVnLX+evBYyUB+t21GfULBpNU4jb/wAwob/R4BWYs8O5cYdviggFkRA9Sgt8Jyiy+PMlh9KjtID5q5K96x+cNjbjPpBsfOPmaHULP3kYDgnPXn1EmmALimVayhSqatG3spia8kIOGtCNMyf246NNe30d1ZB8F6xXItPEWUF47brpnuZCz7o6tLi22uUzvdJdXVhogtPpckMppfYAT2FENPUBzgTeFpGvZ9nveaXUSKXUMKXUT8x1N5uCAmXwbaXUWKXUgUqpx3KpvDHOwrtRTsk667iJlvPTi9vOPNB3m9eDcPaqLT9IQSRkpEN3CDG3X8GKt7fyHTnpXlrAbWceyJ8vOTRt29zrpnvWrdAUFtk4amRPDuhXaZt5rLk+suFs2JypQTI1swf2r2z1HAnWwLwgPovupQUsuGkmEwZ0o6IomjbK3o1X6GxbGFxdapu4sh3T+R609cP+wWnj6F1RRO/K9A6HNuN8fshkQWlvgvgsTgcuBoYDfwYOVUptEZESYBnw67zWMGPd/E1K4RTNIrneT/2fNryKKUP9Gxovc4tX77AgHDIz3DpmdTt8MCcd2Icv3PcmAN8+fiQje5UxfVTPtP0BLpgy0DPz7ZBqb0FXHA1T1QozVFDn+JCqUi49YggXTBlIZXGUS2Yv5L31u20t4y+XHsqehmD+kiBYmkEQn0WuWO9Fr4oi1u9s4FvHpQ8qbC25mE9aa2oZXF3Cpj2NnHhAb2ZNHZTlHK06haYL4XyPOsN8Fl8EfqGUmudcqZSqF5FL81Ot4ATyWQSIRQ6JZLzZmcxQTgoiIeKuPFTRcIhBpmN6QI9iygojaWkr3OSSDuJnZ49nQPcSe3a7bOQqLEIh4funJs1VVx09jKv++rZtwjlyhLfQays9K9o+z4Mbq69QHA1Te/sp7XpsKzw5n9Nb/m7WJBbV7kwLWd4fGN27gtrt9SnT72oyk6pZdHzo7C2APUuPiBQDvZRStUqpV/JVsaAEiYZyquV+ufXDIcn4kXuaoTweTmEkZEzh6hAWkbDQvbSA750yxp4LIRu5qJTHjKrJXshBQThs/m/dy2XVra0mHD9WbjHyWk0yHeTtST5V9a8dM4yWeMKeH92L2RcfwuuO+c1zpVtJQdb5sq1L7GqN7l3nHsRFGwZ7mtY03nQ2zeJJwDkRc9xcd0heapQj/jPlJX87BYefQAhn0yw8BIOX4LE0C6cZytISLjtyqO/x/RjXt4KlG/fkvF8mctUs3FjX04qgskCs2W6MyZg0qP2FheQxvV5JQYQbTx6Tscz0UTVMz1G450rP8kKuO34kp7nm7ujslBZGmDq0qqOr0WUJ53kwYxBhETHTdQCglGo2Q2E7HEECzcEdJH9K2CNvkBOv/bzKF4RDKJUcLZ7pnNl49dqj6VleyIG3zGnV/n7YwqLNmkV+pMVDFx/K6x9tbZfpRt18Huz4IsLVx7afL0bTubHG3OQ7p1YQYbFVRE63IphE5Ayg9Xp0O+M7U56Pn8LPFxAOZdYsvGKYvY5lNcTOUeKttWEP9YiWykaQEdmWkGitZmH5b9oadurH0SN7cvTI/PhBNJr9DREBpTo+Ggq4EviriNyLEa26DvhKXmuVA37jLJxCNsWu56OqeWUkzYZnNJTZADtn28t1DuLW8uq1RwcKV7Wus9XCwjTJ5ctnodFoghMSwzeQ75DpIIPyPgGmikiZubwvyy6fHZIpGsrPDJXJZ9H26R8t53FTihmqbVEK/+/EUYFsuUE1EctE1lYzVL40uuynUAAAHHBJREFUC41GExzDD6fIYxAeEDCRoIicAowDiqxeulLqh3msV2CCZJ0NB4iG8poYJxteQsDWLBxTs7Z1tO7Xpg/PXigHWhwDCFuDdQ+7oqiw5FtXrLtG44XVvHW4ZiEivwNKgGOAB4GzgQV5rVVAjBHcPoPyUhIJZo+GCmXxWXiR0QwV++zNUEE5fHgVZx3cj28fP7JV++c7dDafWD6dcX3T54rWaLoitrDoBD6Lw5VS40XkPaXUrSJyF/BCXmuVA35zJvhNCuIbDSWS8832G2cBqT6LfOeZz5XCSJi7vzSh1fsnQ2e7nrTo262Yp686jLF99HzOmv0DS6PId0qwIMLCSn9aLyJ9ge0Y+aE6Bbn7LPwd3Lnil+4DUs1QnU2zaCtdWbMAmDQo96k+NZrOitW65Dt0NkiX91kR6QbcAbwN1AKP5LNSueA1fSm4Zspzhs667PRTh/agW0mUi6cNzvncXnmmvMxQnU2zaCsR7eDWaDod+U7tnlGzEJEQ8Io5c93TIvJPoEgptTuvtQqIZIqG8kv34dIgBvYo4bErDmvV+b1GTHoKi/1Ms7DngtCyQqPpNHToHNzm7Hj3OZabOougsAgygjt1gF7qJfcobf0oYbdJKxwS+1zPvrvRXp9LUsCuQEhrFhpNp8E5e2I+CdKKvSIiX5ROOn2V3wjuFDOUM0W5o5d//qED+eYM77QIVx49LOu53UIgHPJ2knf1aS7dWJejZYVG0/HYsyd2gsmPvoqROLBJRPaIyF4Rad/Mdq0kcG6olESCyd/XzBxBsc8EO+P7Z4+WcT8bY6xGerl8pqzuCKxkfFqz0Gg6Hmuy43yPs8jaiimlypVSIaVUgVKqwlzusCD1biWp6SwCZZ31NUn539zjxvbiosMHZ6yLW9kKh8Tzge1v01xaz+DiaUM6uCYajcbqs3V4IkEROcprvXsypM8Kd4ppv9xQTpwde/ERHG6i4RC3nD6O2f+pDVy3SMg7Zcj+JSqMmezae+IgjUbTOqwWMN+OgiDjLL7j+F0EHAosBo7NS42y4L4hfpqF00KSMh+3j3nKojAS4vBhrcupHw6F8p6fRaPRaJxYUwV0eNZZpdRpzmURGQDck7ca5UCmObid9vRcBuit+PFJra5PxMcMpdFoNPnCdnB3tM/Cg/VA5unAPkP8oqGcMsRpbnLez/aWxNkmUNJoNJr2xg6d7WjNQkR+7awPMAFjJHeHIK7ffpqFcxY3v7kt2rthj4S9Q2f1nMIajSZfWG1dZ8gNtcjxOwY8qpR6M0/1yRm3z6JHaQE76ppTUlA7hYJfZJQfPz97PD0DTu/pHJQHxujwZ79+RKAJiTQajaY1WG1dh0dDAU8BjUqpOICIhEWkRClVn9eaBcSpWdx1zkH8d/V2nli0PkWz8PNZBFHbzp08IHBd3HNizPt/xwTeV6PRaFpDZ/JZvAIUO5aLgZfzU53cUKRqFtOGVzsGjCXLpYbO5q8+4VBIO7g1Gk2H0BlGcBc5p1I1f5fkr0qZcfoc4gmVIixCAtYwhyChs+1Na2bb02g0mvagM+SGqhORidaCiEwCGvJXpeBY04NaGEIhPRVFajRU/u5o2Cfdh0aj0eSbDh9nAVwDPCkiGzFa4t7Al/Jaq4C4A6FC4khyl7Le22fR3uhxFhqNpqPId67XILmhFgKjgauAK4ExSqnFea1VBjLdjpCI7ZNIcXA7fRbtXJ93bz7eziEV0mYojUbTQXT4THki8r9AqVLqA6XUB0CZiHwtr7VqJSERz4l5vFJ8tNdUp5UlUbqXFABGqhDrXPtbWnKNRtO56Qw+i8vNmfIAUErtBC7PX5WykOGGSCgpGFJ8FqF0n0V7agAxcxT5mD4V9nG1hqHRaD5LOkM0VNg58ZGIhIGC/FWp9YREOGdyfwBmjO5lr09NV278b895sZdv2gvAuL5JYaE1C41G81nSGcZZvAg8LiIzRGQG8CjwQl5rlYHMPgsY17eS2ttPYWBViWO9Q7Og/Xv+Y/sY03scOqSHfa58S3mNRqNxku9IzCDRUNcDV2A4twHew4iI6nT4SdawRzRUe/ksAL5+7HDOO3QAfSqLWb2tDtCahUaj+WzpDNFQCeC/QC3GXBbHAh8GObiInCgiK0RkpYjckKHcF0VEicjkYNX2O473+pDHoLz21Cwi4RB9Ko1B7mH7+HrAhUaj+ezosNxQIjISON/82wY8DqCUCpTwyPRt3Acch5HWfKGI/EMptcxVrhz4JoZAahN+mkXII91He/osvM6lNQuNRvNZ0pE+i+UYWsSpSqkjlFK/BuI5HPtQYKVSapVSqhl4DDjDo9yPgJ8BjUEOmknV8hUWHiO48xWtZM+Hq4WFRqP5DMm3MSPT4c8CPgXmisgDpnM7lxawH7DOsbzeXGdjphEZoJR6LtOBROQKEVkkIovq6up8y/m1z6kpyo3/kXb0WTixclXl6/gajUbjRYdpFkqpZ5RS52GM3p6LkfajRkR+KyLHt/XEIhIC7gauzVZWKXW/UmqyUmpyaWlppmP6rE/+zvegOev4NeXB5sDQaDSa9iDf1owgDu46pdQj5lzc/YElGBFS2dgAOCeD6G+usygHDgBeE5FaYCrwj7Y6ub3wmlY1Xw7ogVUl3Hbmgfx21qS8HF+j0Wi86Awz5dmYo7fvN/+ysRAYISJDMITEecAFjmPtBqqtZRF5DbhOKbWIDLTmfnhFQ+XTAX3BlIF5O7ZGo9F40RkG5bUKpVQMuBp4CSPU9gml1FIR+aGInJ6v83oRCnlpFtqnoNFo9h/yLSxy0ixyRSn1PPC8a93NPmWnBzpoG++HdUPbc1CeRqPRdDQd7rPobLT1doS0ZqHRaPZD8j2VTpcTFm1FbJ/F5+7SNRrNfozWLNoZ63ZqzUKj0exPdFkHd75or9uhfRYajWZ/QguLNNp2Q6xJkbRmodFo9ic6w0x5nYs23pCYlY5D+yw0Gs1+hPZZ+FAYaV3VY3GtWWg0mv2PLj3OIh8IoICSgjBNsUTg/X5+9njG9K5gxWZjClSdQlyj0exPdIY5uDsVZgZwKoqjOe137uQBHNi/knjCEDBas9BoNPsT2mfhQpkO6sochYVFTKcQ12g0+yH5nimvCwoL43+rhYX2WWg0mv2QDp+Du7NhCYtczVAWOhpKo9Hsj+hoKBcJ0s1QpxzYh+qyYJMNWT4L7eDWaDT7E51qPovOgJcZ6lfnHxxYqk4ZUgXAsaNr2r1uGo1G01HkOxqqCwoLQ1pUFCWFRS736KAB3Vh128l5v7EajUbzWaLTfbiwQmfLipJyLlfHjhYUGo1mf0NHQ7mwzFCF4S5XdY1Go8kbej4LF5YZqqCV6T40Go1mf0RHQ7mwQmYP7F/ZwTXRaDSazoPODeWisjjK/J+cRESboTQajcZGp/vwQAsKjUajSUWP4NZoNBpNh6OFhUaj0WiyooWFRqPRaLLS5RzcTh65bAqf7m7s6GpoNBrNfk+XFhaHD6/u6CpoNBrN5wJthtJoNBpNVrSw0Gg0Gk1WtLDQaDQaTVa0sNBoNJouTHnRZ+N67tIObo1Go/m88+q109le15T382hhodFoNF2YnuWF9CwPNq10W9BmKI1Go9FkRQsLjUaj0WQlr8JCRE4UkRUislJEbvDY/m0RWSYi74nIKyIyKJ/10Wg0Gk3ryJuwEJEwcB9wEjAWOF9ExrqKLQEmK6XGA08BP89XfTQajUbTevKpWRwKrFRKrVJKNQOPAWc4Cyil5iql6s3F+UD/PNZHo9FoNK0kn8KiH7DOsbzeXOfHpcALeayPRqPRaFpJpwidFZFZwGTgaJ/tVwBXAAwcOPAzrJlGo9FoIL+axQZggGO5v7kuBRGZCdwEnK6U8hxZopS6Xyk1WSk1uWfPnnmprEaj0Wj8yaewWAiMEJEhIlIAnAf8w1lARA4Gfo8hKLbksS4ajUajaQN5ExZKqRhwNfAS8CHwhFJqqYj8UERON4vdAZQBT4rIOyLyD5/DaTQajaYDyavPQin1PPC8a93Njt8z83l+jUaj0bQPncLB3VZaWlpYv349jY16ilWNP0VFRfTv359oNNrRVdFouhz7hbBYv3495eXlDB48GBHp6OpoOiFKKbZv38769esZMmRIR1dHo+ly7Be5oRobG6mqqtKCQuOLiFBVVaW1T42mlewXwgLQgkKTFf2OaDStZ78RFhqNRqPJH1pYtBPhcJgJEyYwbtw4DjroIO666y4SiQQAixYt4hvf+AYATU1NzJw5kwkTJvD444/z+uuvM27cOCZMmEBDQ0Pe6jd79mw2btxoL1922WUsW7Ys4z7Tp09n0aJFAAwePJht27a1e71uvvlmXn755Vbte88991BfX5+9oEajaTP7hYO7M1BcXMw777wDwJYtW7jgggvYs2cPt956K5MnT2by5MkALFmyBMAue+WVV3LjjTcya9asQOdRSqGUIhTKTc7Pnj2bAw44gL59+wLw4IMP5rR/PojH4/zwhz9s9f733HMPs2bNoqSkpB1rpdFovNCaRR6oqanh/vvv595770UpxWuvvcapp57Kli1bmDVrFgsXLmTChAn8/ve/54knnuD73/8+F154IQB33HEHhxxyCOPHj+cHP/gBALW1tYwaNYqvfOUrHHDAAaxbt8633JgxY7j88ssZN24cxx9/PA0NDTz11FMsWrSICy+80NZgnFrDVVddxeTJkxk3bpx9LD9uvvlm7rnnHnv5pptu4pe//GVKmdraWkaPHs2FF17ImDFjOPvss20NYPDgwVx//fVMnDiRJ598kosuuoinnnqKF198kXPOOcc+hnXP/Or3q1/9io0bN3LMMcdwzDHHADBnzhwOO+wwJk6cyDnnnMO+ffta9wA1Gk0a+51mceuzS1m2cU+7HnNs3wp+cNq4nPYZOnQo8XicLVuSWUxqamp48MEHufPOO/nnP/8JwFtvvcWpp57K2WefzZw5c/j4449ZsGABSilOP/105s2bx8CBA/n444/505/+xNSpU7OWe/TRR3nggQc499xzefrpp5k1axb33nsvd955p63hOPnJT35Cjx49iMfjzJgxg/fee4/x48d7Xtcll1zCWWedxTXXXEMikeCxxx5jwYIFaeVWrFjBH/7wB6ZNm8Yll1zCb37zG6677joAqqqqePvttwF48cUXAZg5cyZXXHEFdXV1lJaW8vjjj3Peeef51u8b3/gGd999N3PnzqW6uppt27bx4x//mJdffpnS0lJ+9rOfcffdd3PzzTen1U2j0eSO1iw6EXPmzGHOnDkcfPDBTJw4keXLl/Pxxx8DMGjQIKZOnZq13JAhQ5gwYQIAkyZNora2Nut5n3jiCSZOnMjBBx/M0qVLM/oyBg8eTFVVFUuWLLHrUFVVlVZuwIABTJs2DYBZs2bxxhtv2Nu+9KUvpZWPRCKceOKJPPvss8RiMZ577jnOOOOMwPWbP38+y5YtY9q0aUyYMIE//elPrFmzJuu1azSaYOx3mkWuGkC+WLVqFeFwmJqaGj788MNA+yiluPHGG/nqV7+asr62tpbS0tJA5QoLC+3lcDic1Wm+evVq7rzzThYuXEj37t256KKLso5FuOyyy5g9ezabNm3ikksu8SzjDlN1Ljuvxcl5553HvffeS48ePZg8eTLl5eWB66eU4rjjjuPRRx/NWHeNRtM6tGaRB7Zu3cqVV17J1VdfnVNs/wknnMAf//hH29a+YcOGFDNWruWclJeXs3fv3rT1e/bsobS0lMrKSjZv3swLL2Sff+rMM8/kxRdfZOHChZxwwgmeZdauXctbb70FwCOPPMIRRxyR9bhHH300b7/9Ng888IBtgspUP+c1TZ36/9u7/9io6zuO4883rVImUVbc/oDqgPBDCtaWNd2oIaEtIkRSiatFcMgPMwJhEJQFBa1DohkzJiT8yKobUAZs/AbZYBnRoBPDhkgIg/Ij0DSO2ahhUiNRfshnf9z1vJZe7679fns/+nokDdznPt/vfe6d7937Pt/v3fvzUz744APOnz8PwJUrVzh37lzUxxSR2KTdzCJRvv76a/Lz87l+/TqZmZlMnTqVZ599Nq59jB07ltOnTzNy5EgAevbsyaZNm8jIyGhXv3DTp09n9uzZ9OjRI/QmDvDAAw9QUFDAfffd1+zUUVtuv/12SkpK6NWrV8THHDJkCGvWrGHmzJnk5uYyZ86cqPvNyMhgwoQJ1NTUsGHDhqjjmzVrFuPGjaNPnz4cPHiQmpoaJk+ezNWrgWVRXnnlFQYPHhz1cUUkOnPOJXoMcSksLHRN3+Jpcvr0aYYOHZqgEXU9N2/eDH2badCgQbfcX19fz4QJEzh58mQCRtc2HSvSVZnZR865W7/hEiOdhpK41NbWMnDgQMrKylpNFCKSnnQaSuKSm5tLXV1dm3369euXlLMKEWk/zSxERCQqJQsREYlKyUJERKJSshARkaiULDxiZixcuDB0+/XXX2fp0qWdOoZPPvmEioqKdm17/Phx9u/f7/GIRCRdKFl4pHv37uzatcuXNR9icePGDfr06cOOHTvatb2ShYi0RcnCI5mZmcyaNYsVK1bccl99fT2lpaXk5eVRVlbGxx9/DAR+VT1//nyKi4sZMGBAxDf6pl9fFxYWMnjw4FDF2pqaGsrLyyktLaWsrIz6+nqGDx8OBMpfnDp1KrSPppLkR44cYeTIkRQUFFBcXMzZs2e5du0aL730Elu3bg0tynTlyhVmzpxJUVERBQUFvPXWW16HTERSSHr+zmL9I9/9f8a+jreHt7Vh7ty55OXlsWjRombt8+bNY9q0aUybNo1169Yxf/589uzZA0BDQwOHDh3izJkzlJeXRzyNVF9fz5EjR7hw4QIlJSWhGkjHjh3jxIkTZGdnN6swO2nSJLZt28bLL79MQ0MDDQ0NFBYW8uWXX/L++++TmZnJ22+/zZIlS9i5cyfLli3j6NGjrF69GoAlS5ZQWlrKunXruHz5MkVFRYwZMyZiEUARSW+aWXjozjvv5KmnnmLlypXN2g8fPsyUKVMAmDp1arNy3RMnTqRbt27k5uby6aefRtx3ZWUl3bp1Y9CgQQwYMIAzZ84A8NBDD5Gdnd1q/6aZyrZt20JJqLGxkccff5zhw4fzzDPPNJt9hDtw4ADLly8nPz+f0aNH880334RmRCLS9aTnzCLSTMCr9jYsWLCAESNGMGPGjJj6h5cUb6rT9cILL7BvX+Cxm5ZfjVTyO9In/b59+9K7d29OnDjB1q1bqa6uBqCqqoqSkhJ2795NfX09o0ePbnV75xw7d+5kyJAhMT0PEUlvmll4LDs7m8rKStauXRtqKy4uZsuWLQBs3ryZUaNGtbmPV199lePHj4cSBcD27du5efMmFy5coK6uLqY38UmTJvHaa6/R2NgYWvmusbGRvn37AoFrHk1aljB/+OGHWbVqVSiBNa0dLiJdk5KFDxYuXNjsW1GrVq1i/fr15OXlsXHjxlvWrI7FvffeS1FREePHj6e6upqsrKyo21RUVLBlyxYqKytDbYsWLWLx4sUUFBRw48aNUHtJSQm1tbWhC9xVVVVcv36dvLw8hg0bRlVVVdxjFpH0oRLlKWD69OmhdbqlY9L9WBGJRCXKRUTEd+l5gTvNhF9bEBFJhLSZWaTa6TTpfDpGRNovLZJFVlYWly5d0puBROSc49KlSzF9MUBEbpUWp6FycnK4ePEin3/+eaKHIkksKyuLnJycRA9DJCWlRbK47bbb6N+/f6KHISKStnw9DWVm48zsrJmdN7PnW7m/u5ltDd7/LzPr5+d4RESkfXxLFmaWAawBxgO5wGQzy23R7WngC+fcQGAF8Fu/xiMiIu3n58yiCDjvnKtzzl0DtgCPtujzKLAh+P8dQJm1LIIkIiIJ5+c1i77Af8JuXwR+EqmPc+6GmTUCvYFmKwiZ2SxgVvDmVTM76cuIv3MX0OjztrH0i9QnnvaWbS1v302LePsgGeLZ1v2xxjNafDsjlpHG4fV2XsezKx+bsfTtjNd6x6qCOud8+QMqgD+E3Z4KrG7R5ySQE3b7AnB3lP0e9WvMYY/xpt/bxtIvUp942lu2tXK7S8SzrftjjWe0+HZGLDsSz3i28zqeXfnY7Eg8k+m17udpqP8C94Tdzgm2tdrHzDIJZMJLPo4pVn/phG1j6RepTzztLds68tzaKxni2db9scYzlvh2hvY+ZjzbeR3PrnxsxtI36V/rvhUSDL75nwPKCCSFD4EpzrlTYX3mAvc752ab2RPAY865ylZ3+N02R10HimFJc4qndxRLbyme3upoPH27ZuEC1yB+CfwdyADWOedOmdkyAtOhvcBaYKOZnQf+BzwRw67f9GvMXZTi6R3F0luKp7c6FM+UK1EuIiKdLy1qQ4mIiL+ULEREJColCxERiSqtkoWZDTWzajPbYWZzEj2eVGdmE83s98H6XWMTPZ5UZmYDzGytme1I9FhSlZndYWYbgsfkk4keT6qL95hMmmRhZuvM7LOWv86OVowwnHPutHNuNlAJPOjneJOdR/Hc45z7BTAbmOTneJOZR7Gsc8497e9IU0+csX0M2BE8Jss7fbApIJ54xntMJk2yAGqAceENkYoRmtn9ZvbXFn8/DG5TDuwD9nfu8JNODR7EM+jF4HZdVQ3exVKaqyHG2BL4YW9TCaFvO3GMqaSG2OMZl6RZz8I5949WSpSHihECmNkW4FHn3G+ACRH2sxfYa2b7gD/5N+Lk5kU8g0UdlwN/c84d83fEycurY1NuFU9sCdSXywGOk1wfdJNGnPGsjWffyR7w1ooR9o3U2cxGm9lKM3sDzSxaE1c8gXnAGKDCzGb7ObAUFO+x2dvMqoECM1vs9+BSXKTY7gJ+Zma/IzGlQVJVq/GM95hMmpmFF5xz7wLvJngYacM5txJYmehxpAPn3CUC136knZxzV4AZiR5Huoj3mEz2mUUsxQgldoqndxRL/yi23vIknsmeLD4EBplZfzO7nUDtqL0JHlMqUzy9o1j6R7H1lifxTJpkYWZ/Bg4DQ8zsopk97Zy7ATQVIzwNbAuvWiuRKZ7eUSz9o9h6y894qpCgiIhElTQzCxERSV5KFiIiEpWShYiIRKVkISIiUSlZiIhIVEoWIiISlZKFiEfMrLyp/LOZLTWzXyV6TCJeSavaUCKJ1FTxONHjEPGDZhYiQWb2czM7YmbHzewNM8sws6/MbIWZnTKzd8zsB8G+882s1sxOBEs+Y2bTzWx1K/vNN7N/BvvuNrPvB9vfNbPfBh/znJmN6txnLBI7JQsRAkvyElgN8EHnXD6BxXWeBO4AjjrnhgHvAb8ObvI8UOCcyyN65c4/As8F+/47bB8Amc65ImBBi3aRpKLTUCIBZcCPgQ8Daz7RA/gMuAlsDfbZRGBNBYATwGYz2wPsibRTM7sL6OWcey/YtAHYHtalaX8fAf06/CxEfKKZhUiAARucc/nBvyHOuaWt9GsqpvYIgaUqRxBIMO394HU1+O+36MObJDElC5GAdwisCNi0lnu2mf2IwGukIthnCnDIzLoB9zjnDgLPAXcBPVvbqXOuEfgi7HrEVAKns0RSij7JiADOuVozexE4EEwG14G5wBWgKHjfZwSua2QAm4KnmAxY6Zy7HDx91ZppQLWZfQ+oQ6u9SQpSiXKRNpjZV865VmcNIl2JTkOJiEhUmlmIiEhUmlmIiEhUShYiIhKVkoWIiESlZCEiIlEpWYiISFRKFiIiEtX/AXESQyZsTEw/AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "import pickle\n",
    "\n",
    "epsilons, baseline, accuracy = pickle.load(open(\"lr_accuracy_500.p\", \"rb\"))\n",
    "\n",
    "plt.semilogx(epsilons, accuracy, label=\"Differentially private\")\n",
    "plt.plot(epsilons, np.ones_like(epsilons) * baseline, dashes=[2,2], label=\"Non-private\")\n",
    "plt.title(\"Differentially private logistic regression accuracy\")\n",
    "plt.xlabel(\"epsilon\")\n",
    "plt.ylabel(\"Accuracy\")\n",
    "plt.ylim(0, 1)\n",
    "plt.xlim(epsilons[0], epsilons[-1])\n",
    "plt.legend(loc=3)\n",
    "plt.show()"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.8"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
